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The isostaticity theory for stress transmission in macro- 
scopic planar particulate assemblies is extended here to non- 
rigid particles. It is shown that, provided that the mean 
coordination number in d dimensions is d -\- 1, macroscopic 
systems can be mapped onto equivalent assemblies of perectly 
rigid particles that support the same stress field. The error in 
the stress field that the compliance introduces for finite sys- 
tems is shown to decay with size as a power law. This leads 
to the conclusion that the isostatic state is not limited to in- 
finitely rigid particles both in two and in three dimensions, 
and paves the way to an application of isostaticity theory to 
more general systems. 
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Much attention has been given lately to particulate 
systems both due to their overwhelming technological 
importance and the fundamental theoretical challenges 
that they pose [1]. In particular the micro- and macro- 
mechanics have focused research activity following exper- 
imental [2] [3] and numerical [4] observations of nonuni- 
form stress fields [5]. Specifically, stresses frequently ap- 
pear to be supported by arch-like regions, termed force 
chains, that cannot be straightforwardly described by 
conventional approaches [6] . It has been recognized that 
to understand this phenomenon it is essential to first un- 
derstand transmission of stresses in 'isostatic' systems [5] . 
Isostatic states are configurations of particles where the 
interparticular contact forces are statically determinate, 
i.e. from the mechanical equilibrium conditions of bal- 
ance of force and torque moments. This means that the 
interparticular forces can be determined without refer- 
ence to compliance and hence to stress-strain relations. 
Isostatic states are characterized by low mean coordina- 
tion numbers per particle which depend on the dimen- 
sionality of the system and on the particles roughness. 
For rough and infinitely rigid particles in d-dimcnsional 
systems [d — 2, 3) this number is Zc = d + 1, for smooth 
infinitely rigid particles of arbitrary shape Zc — d{d + 1) 
[7] [8] [9] , and for smooth infinitely rigid spheres Zc = 2d. 
Isostatic packings of particles are marginally rigid and 
such states have been shown to be easy to approach ex- 
perimentally [9] , making them interesting more than only 
theoretically. Several empirical [10] [11] and statistical 
[5] [8] models have been proposed for the macroscopic 
stress field equations in these systems, suggesting a lin- 
ear coupling between the components of the stress tensor. 
This has been recently established from first principles in 
the two-dimensional case for systems of infininitely rigid 
particles [12]. The new isostaticity theory (IT) closes the 
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FIG. 1. The loading on a packing of grains by a force Fy 
that is distributed evenly on the shown surfaces 

stress field equations with a constitutive relation between 
the stress tensor a and a rank-two symmetric fabric ten- 
sor P which characterizes the local microstructure: 

Pxx^yy ~^ Pyy^xx 2pxy^xy — ■ (1) 

On the scale of a few particles, this equation is a lo- 
cal manifestation of the torque balance condition beyond 
the global requirement that a — [12]. It then tran- 
spired that the coarse-graining of eq. (1) is not trivial, 
but this was eventually resolved, making it applicable for 
macroscopic systems, albeit with a subtle difference in 
the interpretation of the constitutive field pij [18]. This 
paved the way to several results, most notably it enabled 
a derivation of the general solution for the stress field in 
two-dimensional isostatic granular packings [13]. The so- 
lution turned out to indeed give rise to force chains and 
arches. This, not only gave a firm theoretical basis that 
explains the experimentally observed force chains, but 
also provided a way to predict the trajectories of individ- 
ual force chains. Using these predictions made it possible 
to test the theory by direct comparison with experimen- 
tal measurements. 

However, much controversy surrounds the validity of 
the new theory. In particular, because it has been de- 
veloped for infinitely rigid particles there remained ques- 
tions concerning its validity to general particulate sys- 
tems, whose rigidity is unavoidably finite. The clarifica- 
tion of this point is a crucial first step towards bridging 
between IT and elasticity theory. A detailed examination 
of this issue, both in two and in three dimensions, is the 
aim of this paper. 

Two dimensions: 

Consider a polydisperse planar packing of N particles 
of arbitrary shape and typical area a [14], confined to 



within a square container of dimensions L"^ ^ aN . All 
the particles are presumed to be made of the same mate- 
rial whose elastic properties are known [15]. The packing 
is loaded by an infintcsimally small external compressive 
force Fy on two opposite boundaries, as shown in figure 
1. The load compresses the particles against one another 
slightly and the contacts between neighbouring particles 
consist of short lines. The line contacts, rather than point 
contacts, between particles constitute the main difference 
between packings of compliant and infinitely rigid parti- 
cles. The criterion for 'smallness' of Fy is that the contact 
lines are smaller that the linear size of the corresponding 
particles [16]. For a system of infinitely rigid particles 
to be statically determinate in two dimensions the mean 
coordination number per particle must be = 3 up to a 
boundary-to-bulk correction term. We wish to determine 
whether stress fields that develop in assemblies of compli- 
ant particles that satisfy this condition are also governed 
by the equations of IT. 

An ideal resolution of the issue would be to establish 
whether all the interparticular forces can be determined, 
at least in principle, from balance conditions alone. If 
this is possible then the system is statically determinate 
and isostaticity theory must apply. The main difference 
between the geometries of infinitely rigid and compliant 
systems is that while in the former the interparticular 
forces act at a point of contact, in the latter they are 
continuously distributed along contact lines. Let us ex- 
amine the contact between two touching particles, g and 
g' . The particles press on one another with a force that 
is distributed along the contact line with density 4'{x), 
where x is a length parameter that varies from to I 
along the line (see figure 2). Due to the arbitrary shapes 
of particles this force density need not be uniform. Both 
a force and a torque moment are transmitted through the 
contact and these are given by 

faa' = f ${x)dx (2) 
Jo 

and 

Mffff' = / ,f{x) X p{x)dx . (3) 
Jo 

Here p{x) is the position vector from the centroid of the 

particle (defined as the mean position vector of the con- 
tact points of the particle g, see figure 2). 

So, do we need to determine the entire force distribu- 
tions along the contact lines? Considering relations (2) 
and (3), the answer is encouragingly no. The torque mo- 
ment can be represented by a single force of magnitude 
f^a acting at a point x = Iq that lies between x = 
and a; = ^ (it is straightforward to see that Iq cannot be 
outside this section), and whose location is determined 
by the relation 

X p{lo) = M39' . (4) 
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FIG. 2. The distribution of forces, 0(x-) along the contact 
line between grains g and g' . The contact line is parameterised 
by < a; < Z from left to right. The mean force is located 
at a distance a; = Zo, found from the first and second moments 
of the force density 

Thus, it seems that, at least in principle, we can reduce 
the problem to find the discrete forces /^^ . In mechan- 
ical equilibrium these interparticular forces balance out 
[17] 

9' 

The stress field can be defined in terms of the force mo- 
ments around the particles 

s% = Y.n'' pf\io) . (6) 

g' 

The torque balance condition for every particle amounts 
to the requirement that Sf, ~ S%. With this definition 
the stress within a given region inside the material is 
the area average of the force moments over the particles 
within the region. 

Expressions (2), (3) and (6) suggest that the stress field 
is determined only by the forces /^^ , rather than by the 
entire distributions of the contact forces. It follows that if 
we knew the locations of the points where the equiv- 
alent interparticular forces act then we could map the 
system of compliant particles onto an equivalent one of 
infinitely rigid and infinitely rough particles that contact 
at these points, as illustrated in figure 3. The equivalent 
system would have the same mean coordination number 
and it would transmit the same interparticular forces. 
Therefore, it would also have the same macroscopic stress 
field. Since the stress field in the equivalent packing is 
governed by IT then this would lead to the important 
conclusion that the original packing of compliant parti- 
cles is also isostatic and is indeed described by IT. 

However, there still remains the issue of the adaptation 
of the formalism to the structure of compliant systems. 
In particular, recall that IT relies on the identification 
of the geometric tensor whose components depend 
directly on the positions of the contact points. Thus, the 
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FIG. 3. The equivalent system of ideally rigid particles 
(dashed linos) touch at the point x = along the contact 
line of the original compliant particles (full lines) 



question that we are faced with is whether it is possible to 
identiiy the points along the contact lines, /fj^ . Here wc 
appear to have a problem. To determine the locations of 
these points requires using relations (3) and (4), which in 
turn require full knowledge of (^[x) at the contact between 
every two particles. But this is tantamount to a solution 
of the interparticular forces in the first place. Does this 
mean that we cannot find the equivalent rigid packing? 
Have we reached a dead end? 

Not necessarily. The conundrum can be resolved as fol- 
lows. Let us introduce a judiciously chosen approximate 
equivalent system for which we can use IT to solve for the 
stress field. The idea is to show that the difference be- 
tween the approximate and the true fields diminishes as 
iV~" (a > 0) when the system size increases and there- 
fore that the approximate solution converges to the true 
solution for macroscopic systems. The equivalent pack- 
ing is generated by choosing the forces to act at the 
centres of the contact lines. This requires only knowledge 
of the structure, not the force distributions. Let us con- 
struct the geometric tensor P for the equivalent system, 
using the definition in [12] and coarse-grain it using the 
procedure in [18]. Together with the boundary data, we 
can now determine the stress field using the solution of 
reference [13]. 

The deviation of this solution from the 'true' stress 
field arises from the error in the position vectors that 
point from the centroids of the grains to the location of 

fP^ (lo), the error in the stress field around particle g is 



p^^ (lo)- Defining the error vectors as 5^^ 
given by 



(7) 



where is the area associated with particle g. 

Now, a continuous stress field representation is only 
useful on scales that contain a good number of particles 
M but where M <^ N. For macroscopic description, 
the system is regarded as a collection of continuous such 



units. The error made in the stress within such a unit 
region of area = J2n^^ ~ is 



(8) 



9,9'er 



The error in the stress field is linear in the 5^^ and its 
magnitude depends directly on the correlations between 
these quantities; 
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(9) 



Let us now assume that in isotropic packings the er- 
ror vectors S^^ are random and uncorrelated. This, of 
course, may not be the case since correlations may arise 
from the history of the dynamics that gave rise to the 
structure, as well as from inhomogeneities in material 
properties and granular characteristics that lead to non- 
linear contact lines. However, in the absence of evidence 
to the contrary it is plausible that under small loading 
these effects are negligible. Then the sum in (8) can be 
regarded as a two-dimensional Markovian random walk 
and it increases as 0{M), up to logarithimic terms. But 
the area also increases as 0{M) and hence the entire ex- 
pression decreases as 1/M. We therefore arrive at the 
conclusion that 
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(10) 



This result encouragingly support the idea that the ap- 
proximate and the true stress fields converge in the 
macroscopic limit. But we are not finished yet. Hav- 
ing partitioned the system into N/M basic units of M 
grains, we now face the acid test of the analysis. We 
need to determine the size of the discrepancy between 
the boundary data and the corresponding data derived 
from the approximate field. Using the same rationale, 
it is assumed that in isotropic; systc;nis the errors in the 
stresses (10) in different basic units are independent. It 
follows that the error in the stress field at the bound- 
ary (which is normalised by the total area) is of order 
0{y/M/N) <C 1. For macroscopic systems this error 
is indeed negligibly small. We have therefore reached 
the desired result; in macroscopic packings of compliant 
grains the approximate and the true stress fields are the 
same and can be obtained by solving for the isostatic 
stress in the equivalent infinitely rigid packing. 

To make the analysis even more quantitative, let the 
particles' Young modulus be E and let us assume that 
their local radius of curvature is typically R = a^/a. 
In isotropic systems it is expected that the value of the 
parameter a would be distributed over particles around 
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Xj^Jix. For monodispcrse circular particles the distribu- 
tion of a is almost a (5-function around this value. A 
sensible choice of M is such that there are many parti- 
cles in a unit region on the boundary that arc pressed 
by the boundary loading. With this choice the fluctua- 
tions of the force on the boundary particles can be disre- 
garded and the mean force per particle in the y-dircction 
is Fy/\/N. Two particles in contact exerting a normal 
force /„ on one another deform slightly and according to 
Hertz theory the line contact between them is 



ttE' 



(11) 



long. In this expression R'^ is an effective radius, l/W^ = 
l/Rs + 1/R3' « 2y^, E' = 2(1 - u^)/E and u is 
Poisson's ratio. Substituting for the compressive force 
/„ gives that the width of a typical line contact is 




1/4 



(12) 



The error in the distance between the middle of the line 
and the true point is at most w/2. Thus, using re- 
lations (7) and (12), and taking into consideration that 
there are on average three contacts per particle, the typ- 
ical error made in the computation of the stress around 
any one particle is bounded by 



1-^41 



< 



3wf„ 3^/2 



2a 



t3/4 



3/4 



(13) 



This calculation gives the precise power a = 3/4 with 
which the error between the stress fields decays with the 
number of particles N and completes the proof. 

Recalling that theories for ideal rigid particles predict 
stresses that propagate nonuniformly along arches [5] [10] 
[13], this explains why such force chains are also observed 
in packings of compliant particles [2]- [4]. Moreover, since 
the trajectories of force chains can be predicted in ideal 
packings [13] then the above suggests that these predic- 
tions can be extended to systems of compliant particles. 
It would be interesting to compare these predictions with 
the actual trajectories observed in realistic systems, such 
as those of reference [2] . 

Three dimensions: 

The discussion of the three-dimensional systems follows 
the same rationale. The systems considered here con- 
sist of compliant and non-slipping particles, slightly com- 
pressed under a low external load. The mean coordina- 
tion number is four per particle and the contacts make 
small two-dimensional surfaces. We assume, for simplic- 
ity, that the particles have homogeneous elastic proper- 
ties, in which case the contact surfaces are planar. It 
is straightforward to lift this assumption and extend the 
results to systems of particles with nonuniform proper- 
ties. The difference between such a system and one of 



infinitely rigid particles is that the intcrparticular forces 
are distributed across the contact surfaces. The contact 
surface between particles g and g' can be described by by 
a position vector fP^ {x,y), where x and y paramctcrise 
the surface of the contact. Terming the interparticular 
force density {x,y), the mean force between the two 
particles is 



4>^^ (x, y)dxdy 



(14) 



where s stands for an area integration across the contact 
surface. The mean torque moment on the surface is 

M99' = [ $99' {x,y)p^<'' {x,y)dxdy . (15) 

J s 

From eqs. (14) and (15) we now extract a position vector 
Pq^ by using 



f99' X^o^ = M^S' . 



(16) 



It is straightforward to verify that if the contact surface 
is planar then pf,^ corresponds to a point on the surface 



fp9 (xo^yo). Non-planar surfaces, which may re- 



sult from non-miiform clastic properties in the particles, 
do not pose a limitation on this analysis as long as the 
deviation from the plane is smaller than the size of either 
of the particles in contact. We can now define the equiv- 
alent ideal system by postulating that its infinitely rigid 
particles make contacts at the points p^^ . 

This done, we face the same conundrum as in the two- 
dimensional case; the determination of the contact points 
of the equivalent system poses the same level of difficulty 
as the original determination of the contact force distri- 
butions. For the resolution of this problem we follow 
the same logic as before. We construct an approximate 
system, for which isostaticity theory can be applied to 
determine the stress field, and then we show that the ap- 
proximate field converges to the true field as the size of 
the system increases. 

The contact points of the approximate system of in- 
finitely rigid particles are postulated to be the centroid 
of the contact surfaces of the true system of compliant 
particles, 



ggg 

rapprox 



Is P^^' y)dxdy 
dxdy 



(17) 



These points are well defined from the geometry. The 
stress field in the approximate ideal packing of rigid par- 
ticles is 



a^. = — V f^^' (poa' ) 

n yg A^-'^ y^approx J ^ 

g' '' 



(18) 



where V9 is the volume associated with particle g. There 
are several ways to define the volume V9 such that 
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ysys _ V^, but the precise definition is not essen- 
tial for tfie present discussion. The error in the stress 
around particle g originates from the deviations of the 
true positions of the effective forces from the approxi- 

ppg' _ 

rapprox rO 



mate positions SfP^ = nss „ pt^s 



6a?- = — V f?s <5pss 
u yg f^J 



(19) 



Consider a region F containing M <C particles. The 
volume of the region is of order = ~ VM, 

where V is the typical particle volume. The error in the 
stress over this region is 



5al 



(20) 



Following a similar analysis as in two dimensions we 
end up with a three-dimensional Markovian random walk 

of the vectors SfP^ . Assuming now that in isotropic sys- 
tems these are uncorrelated (at least above some length- 
scale), we conclude that (Scrfj) increases as 0(M^/^) 
and therefore that the error decreases as 0{M^^/^). By 
considering many such regions we are then led to the 
conclusion that the error at the boundary between the 
true and approximate stress fields decreases again as 
Oi^/M/N) ~ 0(L-3/2) ^ ^ 

We can make the analysis again more quantitative by 
taking into consideration the elasticity of the particles, 
their typical size and assuming Hertzian interaction. Fol- 
lowing the line of reasoning as that leading to eq. (13), 
we find that according to Hertz theory the diameter of 
the area of contact between two grains is 



w = 2 



K 



1/3 



(21) 



where /„ is the force pressing them together, 1/K — 
(3/4)[(l - ij?)/E + (1 - At'^)/-E'] is an effective elastic 
constant, and = + 1/R^ is an effective curva- 
ture. Taking a cube of particles and pressing on one of 
its surfaces by a force gives that on average per grain 
there is a normal force of order F^/N^^^. Recalling that 
there are on average four contacts per particles, the error 
in the stress around a particle can be bounded by 



IK-I< 



2wU 
V 



n8/3 



i^l/3 



^4/3^-8/9 _ (22) 



Thus we have demonstrated that, just as in two di- 
mensions, the isostatic solution for the equivalent system 
of ideally rigid particles converges to the true stress field 
of the packing of compliant particles as the size of the 
system increases. 

To conclude, it has been shown in this note that iso- 
staticity theory is not limited to packings of infinitely 
particles. Rather, this theory acn be used to describe 



stress fields in macroscopic systems of compliant parti- 
cles. The only condition that such systems must satsify 
is the same as the one for rigid-particle isostatic systems; 
that the number of contacts per particle \s Zc = d + 1 vn 
d = 2,3 dimensions. 

The next step towards bridging between isostatic- 
ity and elasticity theories involves considering packings 
where the mean coordination number is larger than Zc- A 
suggestion in this direction has been made by this author 
[13] and a detailed formulation of an isoelasticity theory 
will be reported elsewhere. 
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